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Abstract: We present both numerical and semi-analytical results on test-particle acceleration in multiple 
parallel shocks. We apply a kinetic Monte Carlo code and an eigenfunction expansion method to calculate 
the distribution functions for electron populations accelerated in subsequent parallel shocks with speeds 
ranging from non- to fully-relativistic. We examine the levels of particle anisotropy at the shocks and 
discuss the implications for AGN and microquasar jets. 



Introduction 

Fermi acceleration at multiple non-relativistic 
shocks has been dealt with in detail by [8] and 
[10]. However, many objects that are likely to con- 
tain multiple shocks are relativistic, e.g., the jets 
of active galactic nuclei (AGNi), [11] the internal 
shock models of microquasars [6] and both the in- 
ternal shock and reverse shock models of gamma 
ray bursts. [4,9] While the external shocks in these 
sources are inevitably strong, internal shocks can 
be weak, [4,6] either due to hydrodynamic pressure 
in parallel shocks or magnetohydrodynamic pres- 
sure in perpendicular shocks. Strong relativistic 
shocks are known to be capable of producing hard 
particle spectra with N(E) oc E~ 2 - 2 via the first- 
order Fermi mechanism. On the other hand the 
first-order mechanism produces much softer spec- 
tra with weak shocks. [7] However, in the case of 
a pre-existing hard power law population of parti- 
cles then either inclusion of the shock drift mech- 
anism [1] or small-angle diffusion [3] leads to an 
amplification of this power law far above that ex- 
pected from purely adiabatic compression. 

In this paper we consider hard power law spectra, 

N{E) oc E- a , 

created at a strong shock by the first-order acceler- 
ation mechanism, which are swept up by a weak 
shock. There is a number of physical scenarios 
in which this is possible; one such possibility is 
in microquasars and the parsec-scale jets of AGNi 



where electrons accelerated at an external shock 
are advected downstream of it towards the follow- 
ing internal shocks. Another scenario is the colli- 
sion of two expanding shells of plasma where the 
weak reverse shock of the first shell sweeps up the 
high energy electrons produced by strong forward 
shock of the second shell. 

In this paper we restrict ourselves to mildly rela- 
tivistic external shocks, applicable to most AGNi 
and microquasars. However, the internal shocks 
of these sources can still have high Lorentz factors 
due to outbursts of the central engine. 

Methods 

The semi-analytic eigenfunction approach is based 
on solving the steady-state particle transport equa- 
tion for the phase space distribution function: 
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where u is the bulk flow speed in the shock 
rest frame, T = (1 — u 2 )~ x / 2 , /i is the pitch- 
angle and D^^ is the pitch-angle diffusion coef- 
ficient. This equation holds separately upstream 
and downstream of the shock, and boundary con- 
ditions are required to find the full solution. One 
condition is that the distributions match at the 
shock front, i.e., /_(p_,/i_,0) = 0), 
where the plus/minus sign denotes quantities up- 
stream/downstream of the shock. We also require 



Particle acceleration by multiple parallel shocks 



Prof 


Vshl 


Vdol 


ri 


V s h2 


Vdo2 


A 


.3000 


.2293 


3.95 


.7290 


.61933 


B 


.7070 


.5965 


3.70 


.9770 


.95731 


C 


.9500 


.9105 


3.25 


.9989 


.99783 


D 


.9800 


.9619 


3.10 


.9997 


.99934 



"¥hi 



Table 1 : The sets used. The speed for the both shocks, 
v s hi and v E h2, and the corresponding downstream flow 
speeds, v^oi and t>do2, and the compression ratio of the 
first shocks, n, are given in the far upstream rest frame 
(observer's frame); all calculations and simulations are 
carried out in the relevant shock rest frames. Compres- 
sion ratio of the second shock is 3 in all cases. 



the distribution to be bounded infinitely far down- 
stream. The distribution must be given far up- 
stream and it is zero when only injection at the 
shock is considered. We can expand / as 

oc 

f(p,/j,,z)= a l {p)Q i {n)cx^(K l z/T) 

i— — oo 

where (Qj(/i),Aj) are eigenfunction-eigenvalue 
pairs satisfying 

If the far upstream distribution is g(p~), then the 
boundary conditions imply 

ff(p_)+5>r(p-)Qr(/o = 52af(p + )Qt(n + ) 

i>0 i<0 

If the far upstream distribution follows a power law 
with index s\, i.e., g{p~) = A pZ Sl , we can find 
solutions with af(p±) — afp± Sl by multiplying 
the matching condition by Q^(u + + fi + ) and inte- 
grating over [i + for j > 0. The same applies also 
if we one considers only monoenergetic injection 
at the shock. In this case, letting 
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Si is found by the condition dct W = as in [7]. 

We also used a kinetic Monte Carlo test-particle 
simulation 1 for comparison. We injected particles 
in the upstream of a one-dimensional step shock, 
and followed them under the guiding-centre ap- 
proximation until they reach a pre-defined escape 
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Figure 1 : Simple model of the flow profile considered 
in this paper. The shocks are far enough apart that parti- 
cles once a particle crossed the second shock it can never 
return to the first shock. 



boundary. This boundary was set sufficiently far 
away in the downstream to make sure the particles 
have reached isotropy in the downstream plasma 
frame. When a particle crossed the boundary, it 
was either "absorbed to the downstream", i.e., re- 
moved from the simulation, or mirrored mimick- 
ing the case of the particle recrossing the shock af- 
ter returning from far downstream. This was done 
with the help of the probability of return (see, e.g., 
[5]). Also particle splitting was used to improve 
statistics: when a particle reached certain energy, 
it was replaced by two "daughter particles", which 
were otherwise identical to their "mother", but had 
only half of the original statistical weight. 

In the simulation particles were scattered at the 
end of each Monte Carlo time-step. The small- 
angle scatterings are elastic in the rest frame of 
the scattering centres, which are taken to be mov- 
ing with the local flow speed, so the energy re- 
mains unchanged over a scattering in the local flow 
frame. Energy losses were omitted in these simula- 
tions to allow for better comparison with the semi- 
analytical results. 

For the first shock the injected particles have a 
small initial energy and random direction. For the 
second shock we inject the particles accelerated in 
the first shock. The particle properties were mea- 
sured in the local plasma frame, so the downstream 
properties for the first shock simply become up- 
stream properties for the second. However, the the 
shocks in each case (see Table 1) were simulated 
separately, so for example particles escaping up- 
stream from the second shock cannot travel back 
to the first shock. 



1. http : // www .iki.fi/joni. tammi/qshock 
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Figure 2: Simulated energy spectra (multiplied by E ) 
downstream of the first shocks for the profiles given in 
Table 1 . Spectra have been shifted vertically to allow for 
better comparison. Solid lines give slopes obtained from 
the eigenfunction method. 
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Figure 3: Distribution of power-law tail particles at 
both D-profile shocks as measured in the upstream rest 
frame. 

Results and Discussion 

As outlined in the introduction the process we are 
interested in the process where a power-law dis- 
tribution, pre-created at a strong shock, is being 
swept and further accelerated by a weak shock. We 
consider four different velocity profiles (see Fig- 
ure 1 and Table 1), for which the compression ra- 
tio of the first shock is calculated from the hydro- 
dynamical jump conditions for a plasma satisfying 
the Jiittner-Synge equation of state. The compres- 
sion ratio of the second shock is 3 in all cases. 

Energy spectra and pitch-angle distributions were 
obtained using both the eigenfunction method and 
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Figure 4: Particle spectra before (dashed line) and after 
the second shock (solid line), as measured in the local 
plasma frame, for the C case. The dotted line shows the 
slope obtained from the eigenfunction method. 

particle simulations, and there is good agreement 
between the two methods. Figure 2 shows the en- 
ergy spectra a few scattering lengths downstream 
of the first shock for all profiles together with the 
semi-analytically obtained spectral slopes. While 
the pitch-angle distribution is isotropic at this dis- 
tance downstream, it can be highly anisotropic at 
the shock front. 2 This can be seen in Figure 3, 
showing the pitch- angle distribution at the both 
shocks for profile D, as measured in the corre- 
sponding upstream rest frame. 

The spectral index doesn't change at the second 
shock crossing due to the fact the the "natural spec- 
tral index" of such a weak shock is greater than 
the index for the injected spectra. The power-law 
is shifted due to the Lorentz transformation across 
the shock, and additionally due to the further accel- 
eration by the first-order process. The compression 
of the plasma at the shock crossing doesn't affect 
the plain energy spectrum, but needs to be taken 
into account when comparing the distribution func- 
tions. In addition to compression and the change 
of frames across the shock, the shock acceleration 
mechanism can lead to significant additional am- 
plification. Example is shown in Fig. 4, where the 
number of the power-law particles at a given en- 
ergy (in the local plasma frame) is increased by the 

2. This poses difficulties for comparing the /i distribu- 
tion in the upstream, as for a relativistic shock (cf. case 
D) one would need to simulate of the order of 10 cross- 
ings to get one crossing with p, = +1 due to the extreme 
anisotropy. For this reason, the data for simulated parti- 
cles does not extend far from /jf = — 1 in this frame. 



Particle acceleration by multiple parallel shocks 




o c l ^ I 

0.1 1.0 10.0 100.0 

Shock proper speed r, V, [c] 



Figure 5: Amplification of spectrum as a function of 
shock proper speed. See text for details. 

factor of 17.4 (case C), and general behaviour of 
the amplification with respect to the shock speed is 
drawn for an example case (spectral index a = 2 
and compression ratio r = 3) in Fig. 5. 

The peak near \i = +1 in the particle pitch-angle 
distribution at the shock measured in the down- 
stream rest frame is due to particles that cross 
the shock and never return. For the first shock 
these particles are not considered as they are not in 
the power-law part of the spectrum; in the second 
shock, however, these are well beyond the thermal 
population and thus may have some signifigance 
e.g. in radiation modelling, as the highest-energy 
particles can radiate their energy away before they 
have had time to isotropise in the downstream of 
the second shock. As expected, Figure 6 shows 
that this proportion increases as the shock becomes 
relativistic. 

The eigenfunction method and test-particle simu- 
lations broadly agree on the level of amplification. 
However, the eigenfunction method has no cut-offs 
in particle energy due to assumptions about energy 
losses. The simulations show that as well as an 
amplification in the magnitude of the power law 
part of the distribution, the power law extends to a 
higher cut-off energy. This will result in a higher 
peak energy in the synchrotron spectra behind the 
shock than what exists ahead of the second shock. 

While this paper restricts itself to weak subsequent 
shocks, future work will use both presented meth- 
ods to examine re-accleration at strong shocks in 
order to produce results comparible to [8] and [10] 
in the relativistic limit. 
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Figure 6: Pitch-angle distribution, as measured in the 
downstream frame, of the particles for cases A (solid 
line), B, (dashed) and D (dotted) for the second shock. 

In addition to energy losses due to, e.g., syn- 
chrotron emission or adiabatic expansion, the cur- 
rent study omits possible effects due to turbu- 
lence transmission (see Tammi, elsewhere in this 
volume), second-order acceleration between the 
shocks [12], as well as injection of new low-energy 
particles into the second shock. 
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